Energy spectra of vortex distributions in two-dimensional quantum turbulence 
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We theoretically explore key concepts of two-dimensional turbulence in a homogeneous compressible super- 
fluid described by a dissipative two-dimensional Gross-Pitaeveskii equation. Such a fluid supports quantized 
vortices that have a size characterized by the healing length f . We show that for the divergence-free portion 
of the superfluid velocity field, the kinetic energy spectrum over wavenumber k may be decomposed into an 
ultraviolet regime (k » ^"') having a universal k^^ scaling arising from the vortex core structure, and an in- 
frared regime (k <K f"') with a spectrum that arises purely from the configuration of the vortices. The Novikov 
power-law distribution of intervortex distances with exponent -1/3 for vortices of the same sign of circulation 
leads to an infrared kinetic energy spectrum with a Kolmogorov k^^'^ power law, consistent with the existence 
of an inertial range. The presence of these k^^ and k^^^^ power laws, together with the constraint of continuity at 
the smallest configurational scale k a allows us to derive a new analytical expression for the Kolmogorov 
constant that we test against a numerical simulation of a forced homogeneous compressible two-dimensional 
superfluid. The numerical simulation corroborates our analysis of the spectral features of the kinetic energy 
distribution, once we introduce the concept of a clustered fraction consisting of the fraction of vortices that 
have the same sign of circulation as their nearest neighboring vortices. Our analysis presents a new approach to 
understanding two-dimensional quantum turbulence and interpreting similarities and differences with classical 
two-dimensional turbulence, and suggests new methods to characterize vortex turbulence in two-dimensional 
quantum fluids via vortex position and circulation measurements. 



I. INTRODUCTION 

Turbulence in three-dimensional (3D) classical fluids is as- 
sociated with a cascade of energy from large length scales de- 
fined by the details of an energy-forcing mechanism, to small 
length scales where viscous damping removes kinetic energy 
from the fluid. This range of length scales, and the range of 
associated wavenumbers k, define the inertial range of energy 
flux 111]. As shown by Kolmogorov in 1941 id, the energy 
cascade corresponds to a kinetic energy spectrum that is pro- 
portional to k^^^^ in the inertial range. Turbulence in a 3D 
fluid is also often associated with the decay of large patches 
of vorticity into ever smaller regions of vorticity ; this Richard- 
son cascade of vorticity provides an important visual picture 
of the fluid dynamics involved in 3D turbulence lUt]. 

Remarkably, two-dimensional (2D) incompressible classi- 
cal fluids exhibit very different turbulent flow characteristics 
due to the existence of an additional inviscid invariant: in the 
absence of forcing and dissipation, the total enstrophy jil] of 
a 2D fluid is conserved in addition to the fluid's kinetic en- 
ergy iH-l^. The fluid dynamics during forced 2D turbulence 
are highly distinctive when compared with 3D flows: rather 
than decaying into smaller patches, jvorticity aggregates into 
larger coherent rotating structures [liot (see for a more 
detailed picture in terms of turbulent stress imposed on small- 
scale vortices). Accompanying these 2D fluid dynamics is an 
inverse energy cascade in which energy moves from a small 
forcing scale to progressively larger length scales, defining an 
inertial range for energy transport with energy flux in a di- 
rection opposite that of 3D turbulence. Eventually, energy 
is transported into flows characterized by length scales that 
are on the order the system size lfl2ll . for which dissipation 
may occur Additionally, there is an enstrophy cascade, in 



which enstrophy is transported from the forcing scale to pro- 
gressively smaller scales |[l3ll . Thus in 2D turbulence, the ki- 
netic energy distribution contains at least these two distinctly 
different spectral regimes. 

Quantum turbulence lITill involves chaotic flow in a super- 
fluid ifTsUTgll and is often associated with a random vortex 
tangle in 3D lITsll . In general, the quantization of circula- 
tion strongly constrains the velocity fields allowed in quan- 
tum turbulence, which must be irrotational everywhere within 
the fluid, yet inertial ranges with k'^^^ spectral dependence 
are still found in 3D quantum turbulence 1I20I1 . In an incom- 
pressible superfluid (such as Hell), the vortex core diameter 
can be neglected for all practical purposes, inspiring the study 
of point-vortex models of superfluid dynamics. Such a model 
was used by Onsager to first predict the aggregation of vor- 
tices within inviscid 2D fluids, and was the context for his 
prediction of the quantization of vortex circulation in a super- 
fluid ll2l1l . Despite the historical importance of this approach 
in stimulating advances in 2D classical turbulence Ii22il . char- 
acteristics of 2DQT remain little known, due in part to the 
difficulty of achieving the necessary 2D confinement for in- 
compressible superfluids. The increasing relevance of 3D 
turbulence concepts to dilute-gas Bose-Einstein condensate 
(BEC) experiments ll23l - l25ll and recent theoretical work on 
2DQT 119, 26-31 have highlighted the need for a treatment 
of turbulence in 2D superfluid systems that incorporates the 
concept of compressibility from the outset. Motivated by re- 
cent experimental demonstrations of the confinement needed 
to study 2DQT in a dilute-gas BECs lf32[ [33I1 our aim in the 
present paper is to present a new approach to solving some of 
the open problems of 2DQT in the context of such a system. 

In a BEC, the vortex core size is non-negligible, and stems 
from the healing length ^, a scale of fundamental importance 
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in BEC dynamics that is typically about two orders of magni- 
tude smaller than the system size js^l- Compressibility also 
allows for a rich array of physical phenomena in these super- 
fluids; in particular, a vortex dipole ll32ll (comprised of two 
vortices of opposite sign of circulation) can recombine, re- 
leasing vortex energy as a burst of acoustic waves. The oppo- 
site process of vortex dipole generation from sound may also 
efficiently occur. 

Recent theoretical studies of decaying quantum turbulence 
in 2D BECs have shown that when the vortex dipole annihi- 
lation process is dominant it sets up a direct cascade of en- 
ergy over the scales associated with the dipole decay, sug- 
gesting that this annihilation mechanism could prohibit an in- 
verse energy cascade from occurring in a compressible super- 
fluid ll29il30ll . Moreover, enstrophy in a quantum fluid is asso- 
ciated with the number of vortex cores; if vortices annihilate, 
then enstrophy may not be conserved, bringing into question 
the existence of energy and enstrophy cascades in 2DQT, and 
the universal nature of 2D turbulence and its correspondence 
wifli 2DQT. 

The general characteristics of 2D quantum turbulence in 
compressible quantum fluids, including the capacity for these 
systems to show an inverse energy cascade, enstrophy con- 
servation, and vortex aggregation have thus remained largely 
unknown. However, a recent study of the formation of vortex 
dipoles during the breakdown of superfluid flow around an ob- 
stacle in a highly oblate BEC experimentally and numerically 
observed aggregation of like-sign vortices into larger-scale co- 
herent structures ||32| . and found time scales over which vor- 
tex number and hence enstrophy may remain constant. The 
vortex clustering effect inhibits the dipole-decay mechanism 
by keeping vortices distant from antivortices (vortices of op- 
posite circulation), and suggests that an inverse cascade might 
be observed under the right conditions of forcing. System dy- 
namics consistent with the existence of an inverse energy cas- 
cade were indeed found in a recent study of forced 2DQT in a 

BEC m. 

In this article we address 2D quantum turbulence in a com- 
pressible quantum fluid from an analytical perspective. We 
determine the kinetic energy spectra of vortex distributions 
in a homogeneous compressible superfluid in a quasi-exact 
manner via an analytical treatment of the physics of the vor- 
tex core. We are thus able to study the properties of vortex 
configurations and their resulting spectra in BEC. We develop 
a technique to sample spatially localized vortex distributions 
with power-law behavior over a well defined scale range. We 
are thus able to identify the conditions for an inertial range 
in fully polarized and neutral systems. A polarized cluster is 
sampled using a specific exponent for the vortex locations rel- 
ative to the cluster center, which is size and scale dependent. 
The specific radial exponent is shown to determine the veloc- 
ity distribution in the classical limit and we thus identify an 
expanding inertial range with a steepening velocity distribu- 
tion. 

Making use of the universality of the spectral region gen- 
erated by the vortex core, we identify an analytical form 
of the Kolmogorov constant that we test against dynamical 
simulations of the damped GPE. The derivation of the Kol- 



mogorov constant occurs for a highly idealized vortex distri- 
bution. Thus the complex flows generated by real forcing re- 
quire that we introduce a new parameter called the clustered 
fraction, and modify our ansatz to account for imperfect clus- 
tering, based on the universality of the Kolmogorov constant. 
The modified ansatz agrees well with the numerical simula- 
tions of grid turbulence, supporting our analytical identifica- 
tion of the Kolmogorov constant. 

II. BACKGROUND 

The starting point for much of BEC theory is the Gross- 
Pitaevskii equation (GPE), which provides a capable descrip- 
tion of trapped Bose-Einstein condensates at zero temperature 
ll34ll . Our model, outlined below, consists of a damped GPE 
(dGPE) description of a finite-temperature BEC which can be 
derived from the stochastic GPE theory ll35l HH . In this sec- 
tion we develop a link between the dGPE and the classical 
Navier-Stokes equation, identifying a quantum viscosity aris- 
ing from the damping. The corresponding Reynolds number is 
defined in direct analogy with classical fluids. We then state 
some key properties of a single quantum vortex, and define 
the decomposition of kinetic energy into its compressible and 
incompressible components. 

A. Damped Gross-Pitaevskii theory 

The damped Gross-Pitaevskii equation of motion for the 
quantum fluid wave function i}/{r, t) has been obtained phe- 
nomenologically jjTtl , within ZNG theory llssll . and via a 
microscopic reservoir theory ll35l l39l lioll . and we will con- 
sider it within the context of the latter framework, for which 
the full equation of motion is the Stochastic Projected Gross- 
Pitaevskii equation (SPGPE). The SPGPE is derived by treat- 
ing all atoms with energies above an appropriately cho- 
sen energy cutoff fcut as thermalized, leading to a grand- 
canonical description of the atoms below tcut- A dimension- 
less temperature-dependent rate y describes Bose-enhanced 
collisions between thermal reservoir atoms and atoms in the 
BEC. Neglecting the noise, we obtain the equation of motion 
for the condensate wave function (in the frame rotating with 
the chemical potential fi) 

dMr, t) 

= (iy - - ^)'A(r, 0- (1) 

For atoms of mass m in an external potential V{r, t), the oper- 
ator £. gives the GPE evolution: 

X<A(r, t) = (-^V2 + Vir, t) + g\ip{r, Ol'j <A(r, t), (2) 

where the interaction parameter is g - Anfi^a/m, for i-wave 
scattering length a. This equation of motion has been used ex- 
tensively in previous studies of vortex dynamics llssi l39l I4TI 
I42I1 and provides a capable description of dynamical BEC phe- 
nomena. In general the damping parameter is small (y <*c 1), 
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and it is typically much smaller than any other rates charac- 
terizing the evolution. 

Defining the Gross-Pitaevskii Hamiltonian 



He 



I'' {2m 



\Vil^(r,t)\^ + V(r,tM{r,t)\' + ^\il,{r,t)\'' 



and condensate atom number 



-! 



(3) 



(4) 



the equation of motion ([T]) evolves the grand-canonical Hamil- 
tonian Kc - He - i-iNc according to 



dKc 
dt 



2y 

n 



r\(jj-£Mr,t)t 



(5) 



The stationary solution minimizing Kc is the ground state sat- 
isfying iJ.il/{){r) = Jjpdir). This is a consequence of the nonlin- 
ear form of the damping in ([T). The damping term arises from 
collisions between high energy atoms that lead to a Bose- 
enhanced growth of the matter wave field, with instantaneous 
energy determined by The equation of motion thus de- 
scribes a system coupled to a thermal reservoir in the chosen 
frame of reference. 

The SPGPE provides a rigorous framework for the dGPE 
derivation, originating from a microscopic treatment of the 
reservoir interaction. In particular, y can be calculated explic- 
itly lf36ll for a system with well-defined reservoir parameters 
/i, T, and e^^u i-e. a system close to thermal equilibrium. In 
essence it is computed via a reduced Boltzmann collision in- 
tegral that accounts for all irreversible s-wave interactions that 
can change the condensate population by interacting with the 
thermal cloud. If the thermal cloud is 3D (i.e. j6 ' = ksT is 
greater than the potential well mode spacing in each spatial 
dimension) the damping takes the explicit form 



r 



where 0[z, s, a\ is the Lerch transcendent, and 



(6) 



(7) 



with AdB = -\/2nfP-/mkBT the thermal deBroglie wavelength. 
The dimensionless rate yo provides a useful estimate of the 
full damping strength when the cutoff tcut is unknown. Equa- 
tion (|6]l is independent of position, and valid over the region 
V(r, t) < 2&jj,/3, provided the potential can be treated semi- 
classically |[36ll . The summation gives Bose-enhancement cor- 
rections due to the Bose-Einstein distributed reservoir atoms, 
and is typically of order 1-20 in SPGPE simulations with a 
consistently determined ener gy c utoff" ll43ll . Typicaly y ~ 

87 



5 X 10 in Rb experiments | 



B. Heuristic derivation of a quantum Reynolds number 

In this section we consider the role of dissipation within 
the dGPE description, and show how to recover the celebrated 



Navier-Stokes equation (NSE). In doing so we find an explicit 
expression for the viscosity which has a microscopic quan- 
tum origin, stemming from s-wave scattering of incoherent 
reservoir particles with a coherent superfluid. While not of- 
fering a practical reformulation (the GPE and its generaliza- 
tions are capable numerical workhorses), this indicates a con- 
nection between the dGPE and the NSE in the hydrodynamic 
regime, allowing the identification of a parameter analogous 
to the kinematic viscosity of classical fluids. 

The fluid dynamics interpretation of the Gross-Pitaevskii 
equation is based on the Madelung transformation, which 
we now apply to the damped GPE ([1]), writing (/((r, f) = 
-y'p(r7f)exp [/0(r, r)], where p(r, f) is the number density of 
the superfluid (number of atoms per unit volume), and 0(r, f) 
is the macroscopic phase of the quantum fluid. The velocity is 
then given by v(r, f) - tiV@{r, t)/m. The resulting equations 
of motion (with implicit t and r dependence) for density and 
velocity are then given by 



|,V.(pv)=^(^-t/e«.), 

ot n 
where an eff'ective potential C/gff is defined as 



f/eff(r,f) 



+ V + gp- 



2m 



(8) 



(9) 



(10) 



The last term is called the quantum pressure, which is very 
small except where p changes sharply, such as near vortex 
cores. By neglecting this term in the absence of dissipation 
we are considering the so-called hydrodynamic regime. 
We now consider the y term in (|9]l 



V ■ Vp 



^ (v(V ■ V) + V 
2 \ p 



(11) 



Note that v ■ Vp s for an isolated quantum vortex. In the ab- 
sence of acoustic energy this will also be a good approxima- 
tion for a system of vortices provided their cores are well sep- 
arated, since the density gradient of each vortex is localized to 
a region where the velocity is dominated by the single-vortex 
velocity field. It should thus be a reasonable approximation to 
neglect the second term in ( fTTT i. In a superfluid the curl term 
in the expansion V(V • v) - V x (V x v) -H V^v may also be 
consistently neglected away from vortex cores; similarly we 
neglect the curl term in V(v ■ v) = 2 (v ■ V)v + 2 v x (V x v) 
when taking the gradient of (fTOl i. We then find that (|9]l reduces 
to a quantum Navier-Stokes equation for the velocity field: 

^+(vV)v = --V(y + ^p) + v,V-v (12) 
ot m 

where the kinematic quantum viscosity is 

hy 

in analogy with classical fluids. In this regime, (fT2l i is coupled 
to the continuity equation 



(13) 



^ + V ■ (pv) = ^(p - Uh) 
ot n 



(14) 
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with hydrodynamic potential 



Unir) = — + V + gp. 



(15) 



The source term in (fl4l drives the system towards particle- 
number equilibrium with the reservoir In the Thomas-Fermi 
regime /i - Unir) ~ 0, restoring approximate particle-number 
conservation. 

Making use of (|7]i, we can give an order-of-magnitude esti- 
mate for the viscosity 



= ^ - '^"^bT 
^ 2m nh 

We can also estimate a quantum Reynolds number as 

UL _ ttH UL 

2a^ ksT 



(16) 



(17) 



for BEC flow with characteristic speed U and length scale L. 
We can write the quantum Reynolds number as 



Re' 



_ dB 



^If, mUL 



fl2 4/7 



(18) 



We note that temperature only enters the expression through 
the deBroglie wavelength of the matter wave field, in the ra- 
tio Ads/a, which is typically very large for a BEC. Note that 
strong scattering corresponds to strong damping, and hence 
low Reynolds number A large deBroglie wavelength corre- 
sponds to a relatively cold system, which is hence expected to 
be weakly damped and have a high Reynolds number. 
We thus have a dimensionless ratio in the form 



Re' 



[deBroglie wavelength]^ [flow momentum] 



[scattering length]- [quantum momentum] 



(19) 



where fi/L is interpreted as the quantum momentum associ- 
ated with the transverse length scale of the flow. 

A concrete example is provided by a recent experimen- 
tal study of 2DQT generated by stirring a highly oblate, 
toroidally confined BEC HJj. The initial system consists of 
~ 2.6 X 10'' atoms of ^^Rb at a temperature of ~ lOOnK. Us- 
ing these numbers in our analysis the dimensionless damping 
parameter yo ~ 6 x 10""* gives a kinematic quantum viscos- 
ity v'' ~ 6 X lO^^yum- s"' . The trapping potential confines the 
flow to an annular channel of width L ~ 30 i-im. The nomi- 
nal flow speed can be estimated from numerical simulations 
of the dGPE Issll . giving a value U - 5/ims ' as the peak 
value occurring in the bulk flow during the stirring sequence . 
These values give an estimate Re^^ ~ 600. We can alterna- 
tively estimate a Reynolds number at the scale of the forcing 
in the experiment, which is of order of the size of the vortex 
dipoles nucleated, d ~ 10^, with healing length ^ ~ 0.5 jum. 
Such dipoles have a characteristic speed v^ ~ 146 /^ms"', for 
which we estimate the Reynolds number of the forcing scale 
as ~ V(/c//v" = 6.2x10"'. These large values suggest that turbu- 
lent flow in a finite-temperature BEC may exist across a wide 
range of length scales if the analogy is made with the classical 
Reynolds numbers that correspond to turbulent flow This 



interpretation is broadly consistent with the experimental and 
numerical observations of chaotic vortex dynamics [Ssll . 

We emphasize that the quantum Reynolds number estimate 
proposed here is applicable to a finite-temperature weakly in- 
teracting superfluid and may provide a general condition in 
analogy with classical fluids that is independent of dimension. 
However, taking the zero-temperature limit gives an infinite 
value for /?e°, and in this regime the superfluid dissipation 
stems from vortex reconnections (or annihilation in two di- 
mensions) and coupling to the sound field ll45l - j47ll . The de- 
tailed description of criteria for superfluidity in the cross-over 
from the zero-temperature to high-temperature regimes is an 
open problem ll48ll . Furthermore, due to the reversed direction 
of energy transfer in 2DQT, the scales of interest have to be 
reexamined; we do not pursue this here. Our aim is to estab- 
lish a conceptual link between the dGPE and the NSE, given 
by Eq. ( fT2l i. In doing so we have shown how to identify the 
equivalent viscosity in a finite temperature BEC. 



C. Two-dimensional vortex wavefunction 

In the remainder of this work we limit our analysis to ho- 
mogeneous compressible quantum fluids in two dimensions, 
and redefine our spatial and velocity coordinates accordingly: 
r = {x,y) = r(cos0, sin0) and v = (v.v, v,,). We thus confine 
our attention to the regime of an eff'ective 2D GPE, with modi- 
fied interaction parameter While 2D BEC systems can be cre- 
ated through extremely tight confinement in one dimension, a 
regime of eff'ective 2D vortex dynamics can also be accessed 
in less oblate systems, giving a 2D analysis wider applicabil- 
ity. For example, although the BECs of References Il32ll and 
VM were three dimensional, the confinement along one di- 
mension was strong enough to limit vortex motion to a plane 
and suppress vortex bending and tilting away from the tight- 
trapping direction. Aspects of BEC dimensionality in regards 
to vortices and Kelvin waves were analyzed in [44], further in- 
dicating that sufficiently oblate 3D BECs may be considered 
2D as far as vortex dynamics and turbulence are concerned. 
At the same time, such systems can remain far enough away 
from the quasi-2D limits in which a Berezhkinski-Kosterlitz- 
Thouless (BKT) transition has been observed Ii49i] , and BKT 
physics may thus be neglected. 

For our analysis of kinetic energy spectra, we require cer- 
tain properties of a quantized vortex, namely the asymp- 
totic character of the wavefunction for large and small length 
scales. The Gross-Pitaevskii equation describing the homo- 
geneous (y = 0) 2D Bose gas is obtained from ([T]i by taking 
y = and using an interaction parameter g2 where / is 

the characteristic thickness of the 3D system OOll : 



diff(r, t) 
dt 



i-^+goMrjfUirJ). (20) 



For example, in a system with harmonic trapping in the z- 
direction characterized by trapping frequency w^, the length 
scale is I — y/^jrh where /, - ^JWfmaK is the z-axis harmonic 
oscillator length, and the confinement is assumed sufficient 
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to put the wavefunction into the z-direction single-particle 
ground state. 

For solutions with chemical potential containing a single 
vortex at the origin (with circulation normal to the plane of the 
quantum fluid) we can write ifsoll 

^i(r,f)= ^J^Qe-^'^^'l''xi'■|^)e^^" (21) 

where ^ = fi/mc is the healing length for speed of sound c - 
-sjfi/m, and «o = /i/gi is the 2D particle density for r » ^ and 
is taken to be a constant. The vortex radial amplitude function 
Xio'), where cr - r/^ is a scaled radial coordinate, is a solution 
of 



(22) 



The boundary conditions are x(Q) - 0, and the derivative 
x' = dx/dcr evaluated at cr = must be chosen such that 
it is consistent with;if(oo) - 1 and;i''(oo) - 0. The value 



A = x'{0) = lim 



^Jn^ 



dr 



is determined numerically to be A - 0.8249 . 
(ISTT i has the velocity field of a quantum vortex 



v(r) = — (+ sin0, 
mr 



: COS 6). 



(23) 
The state 

(24) 



D. Kinetic energy decomposition 

We make use of the decomposition of the kinetic energy 
into compressible and incompressible parts ll28[|5ll . The 2D 
case of the Gross-Pitaevskii energy functional Q can be de- 
composed as £ = Ef^ + Ey + El + Eq, where 



Ek 
Ev 
E, 

E. 



Q - 



- J d^rp{r,t)\y{r,t)\\ 
J dh p(r,f)V(r,f), 
f / d\piv,tf, 
^ J d\ \Vyf^)\\ 



(25) 
(26) 
(27) 
(28) 



Respectively, these define the components of energy that 
can be attributed to kinetic energy, potential energy, inter- 
action energy, and quantum pressure. We are interested in 
the kinetic energy, Ek- We define a density-weighted ve- 
locity field u(r, f) = -sjpir, t)\(r, t), then decompose this into 
u(r, t) - u'(r, f) + u'^(r, f), where the incompressible field u' 
satisfies V ■ u' =0, and the compressible field u'^ satisfies 
V X u*^ = 0. We can further decompose the kinetic energy as 
Ek - Ei + Ec, where the portion of Ek attributed to compress- 
ible or incompressible kinetic energy is defined as 



Ecj = f J^'r lu''''(r,Ol'. 



(29) 



The compressible component is attributed to the kinetic en- 
ergy contained in the sound field, while the incompressible 



part gives the contribution from quantum vortices. Our analy- 
sis below only involves Because we focus on vortex con- 
figurations at instants in time, we drop the explicit time de- 
pendence from the remainder of our expressions. 

In fc-space, the total incompressible kinetic energy £, is 
given by 



where 



r/k) = i- Jt/^r e-'" 'u}(r). (31) 

The one-dimensional spectral density in ^-space is given in 
polar coordinates by integrating over the azimuthal angle to 
give 



mk 



Ei(k)^--^y f d<Pk\rj{k)\^ 

2 i^v^O 



(32) 



which, when integrated over all k, gives the total incompress- 
ible kinetic energy £, - dk Ei{k). 



m. INCOMPRESSIBLE KINETIC ENERGY SPECTRA 

A. Incompressible kinetic energy spectrum of a vortex 

We now consider the kinetic energy spectrum of a single 
quantum vortex in a 2D BEC. For an arbitrary wavefunction 
the decomposition into compressible and incompressible parts 
must be performed prior to carrying out the transformation to 
the spectral representation. However, for a quantum state con- 
taining a single vortex and no acoustic energy [i.e. the single 
vortex wavefunction ( l2Tl i1 we note that the wavefunction 
is automatically incompressible, i.e. the compressible part is 
identically zero: 



y/p(r)v(r)) = V ■ 



/p(r)V ■ V s 0. (33) 



The first term vanishes due to the orthogonality of the density 
gradient and velocity of a vortex, and the second due to the 
form of ( |24] |. Thus the incompressible spectrum is the entire 
spectrum for a single quantum vortex. 

For a single vortex we can thus ignore the incompress- 
ible decomposition and cast the kinetic energy spectrum in 
terms of the properties of the radial amplitude function ;if(cr) = 
-y'p(cr^)/«o obtained from (|22] |. We have 



ih d r"" 

m dk Jo 



d'-r e 



/kr 



■ sin Q 



dr Jo{kr) 



m k^ 



f 

Jo 



da x'icr)Joik^(r), (34) 
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where Jq is the zeroth-order Bessel function of the first kind. 
Similar analysis gives 'Fyik) = -'F^{k). We can thus find the 
one-vortex spectrum [see (|32] |1 

Ej(k)^Q^'F(kO, (35) 
where we define the dimensionless core spectral function 

and we have introduced the unit of enstrophy 

IntP'nQ 



(36) 



n = 



(37) 



giving as the natural unit for the kinetic energy density. 
The core spectral function has the small-A;f asymptotic form 



ki:< 



=1 [Jo 



(38) 



For » 1, Jo{k^o-) is highly oscillatory except at cr = 
where it is unity, and the Taylor expansion of x'icr) can be 
truncated at zeroth order to give 



F{kO: 



/tf»i 



A2 



do- Jo{k^o-) - 



A' 
(^ 



We thus have the asymptotic spectra for a single vortex 



k^ 

,2 



Ejik)\ ^A"- 



(39) 

(40) 
(41) 



The <s; 1 regime arises purely from the irrotational veloc- 
ity field of a quantum vortex, while the k^ '» I regime is a 
property of the core of a compressible quantum vortex. The 
fc^ » 1 regime explicitly depends on the slope of the wave- 
function at the core of a vortex. The cross-over between these 
regions occurs in the vicinity of A:^ a; 1, hence we take k^ - \ 
as distinguishing the infrared {k^ < 1) and ultraviolet {k^ > 1) 
regimes in the remainder of our analysis. The scale k^ - I 
thus serves to define an important length scale of the problem, 
namely = 2n^. In Fig. [T] we see that at this distance from 
the vortex core the deviation of the amplitude from the back- 
ground value is very small. This is the scale beyond which the 
details of the core structure are no longer important in charac- 
terizing the wavefunction, or equivalently that the fluid den- 
sity has approximately reached its bulk value. The irrotational 
velocity field in Eq. (l24l i is the only remaining signature of a 
vortex at this range from its center and beyond. We note that 
our derivation of the k^^ power-law stemming from the quan- 
tum vortex core structure is consistent with recent analysis of 
the Kelvin-wave cascade in 3D Is^l . 

Now that we have identified the properties of a single vor- 
tex, it is natural to ask whether a unit of enstrophy can be 
attributed to a single quantum vortex, and to compare this 
with the quantity defined in (|37T i. The point-vortex model sug- 
gests that this can be done, but gives a singular result, which 




FIG. 1 . Amplitude of the wavefunction for a single vortex solution of 
the Gross-Pitaevskii equation. The numerical solution of ( I22t (solid 
line) is compared with the ansatz ( I44t (dashed line). The inclined 
dashed line shows the slope A of the exact solution at the origin. The 
vertical line is the point a = In. 



is nevertheless known to be proportional to the number of vor- 
tices ll30ll . The problem is also evident if we attempt to eval- 
uate the enstrophy of a single vortex from the spectrum (l35t . 
Multiplying by k~ to produce an enstrophy spectral density 
and integrating this over the ultraviolet regime k^ > 1, we are 

j-OO 
J dklk. 

In this work we therefore define a new asymptotic quantity 
with units of enstrophy as 



^ = lim k'Ei(k). 

k^oo 



(42) 



This quantity plays a fundamental role for a compressible 
superfluid because it completely specifies the large-fc region 
of the incompressible kinetic energy spectrum. Because the 
spectrum in this region of A:-space is determined by the core 
structure of quantized vortices, we call this unit the onstrophy 
to both recall Onsager's contribution to our understanding of 
quantized vorticity in a superfluid and emphasize the differ- 
ence between enstrophy in classical and quantum fluids. For 
a single quantum vortex we find, using ( [35T l, the onstrophy 



^1 



which differs from by the factor A^ - 0.6805 . . 
erty of the vortex core in a compressible superfluid. 



(43) 



a prop- 



B. Vortex wavefunction ansatz and kinetic energy spectrum 

To study 2D kinetic energy spectra we will make extensive 
use of an algebraic ansatz for the wavefunction of a single 
vortex in a homogeneous superfluid. Numerical evaluation 
of the exact core function ( [36] l is not straightforward due to 
the highly oscillatory integrand, and the need to determine the 
vortex amplitude ;t'(cr) extremely accurately over a large range 
of length scales. In order to accurately represent the spectrum 
it will be crucial that our ansatz have the correct asymptotic 
properties for small and large length scales described imme- 
diately above (l23T l. Making use of the slope at the origin com- 
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puted for the exact solution in 
function: 



we use the ansatz wave- 



(44) 



+ (A-i^)2 



The general form of this ansatz has been previously used to 
describe the shape of a vortex core llsoll . but here we use 
a length scale A"'^ that enforces matching the slope of the 
ansatz density distribution to the exact value at the center of 
the core. The state ( l44l i has the irrotational velocity field of a 
quantum vortex specified in jTM and reproduces the asymp- 
totic slope of the exact solution near the origin, as shown in 
Figure [1] 

We now compute the kinetic energy spectrum for a single 
vortex by evaluating ( ISTT i using the form fl4t . Taking A"'^ = 
b for brevity, we have 



in V«o c/ r°° dr Joikr) 



dk 



f 

Jo 



- i- 



h -^)b 
2m 



.,fk(f)-.(f).,(f 



(45) 



where Ij and Kj are modified Bessel functions of the first and 
second kind, respectively, of order j. Since \'Fy?'-\Tx?', we 



find the incompressible energy spectrum of a single vortex 



E\(k) = Q^VA(fc^), 



where 



and where we define 



(46) 



(47) 



f(z) = (z/4m(z/2)Ko(z/2) - h{zl2)K,(zl2)f. (48) 
The function /(z) has the following asymptotics: for z <sc 1 

/(^) = 7+(r + ln(|))f + .-., (49) 



where j 
z »1 



0.57721... is the Euler-Masceroni constant; for 



f{z) 



- + - + . 

z^ z^ 



The function F^ikO thus has the asymtotics 



1 

A^ 



(50) 

(51) 
(52) 



which are identical to those of F(k^). The two functions are 
very similar, with only small differences evident in the cross- 
over region A:^ ~ 1, as seen in Figure|2] We use /^a(A:^) instead 
of F{k^) for describing the kinetic energy spectrum for a vor- 
tex core in the remainder of this work as it is numerically ex- 
pedient and does not alter any of the physical consequences of 
our analysis. Towards the end of this paper we will compare 
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10 



10 



FIG. 2. Comparison of the numerically computed kinetic energy 
spectrum [(135), blue solid line] and that obtained from the core func- 
tion F\(k^) r(l46t, red dashed line] in the cross-over regime ~ \. 
The asymptotic expressions ( 138) and ( 139) are shown by the black and 
green lines respectively. 



the asymptotic results of our analysis with spectra determined 
numerically from forced dGPE dynamics. 

The spectrum of a single vortex is shown in Figure |3] and 
compared with the spectrum of a vortex-antivortex pair (a vor- 
tex dipole), and that for two vortices of the same sign (a vortex 
pair). These two- vortex spectra are analyzed in the following 
section. 



C. Two-vortex spectra 



Extension of the discussion in Section lTlI Al leads us to con- 
clude that a wavefunction that only contains vortices (i.e., no 
sound field) separated by more than a few healing lengths 
will thus be approximately incompressible according to the 
decomposition. The approximation breaks down through the 
non-orthogonality of v and V ^jp(r) near a vortex core due to 
the velocity field induced by the other vortices. However, in 
the close vicinity of a vortex core, where V -s/pir) is signif- 
icant, the velocity is dominated by the velocity field of that 
vortex core. An arrangement of vortices separated by more 
than a few healing lengths will thus be approximately incom- 
pressible. In the following analytical treatment we will ne- 
glect any compressible part that arises from an assembly of 
vortices described by the ansatz (l44l i. 

A two-vortex state in a homogeneous system with no 
boundaries has density-weighted velocity field -y/p(r)[vi(r) + 
V2(r)], where \j is the velocity field around vortex j - 1,2 
taken separately. If the vortex cores are separated by li » ^ 
it will also be a very good approximation to write p(r) = 
Pi(r)p2(r)/no, with 



no |r ■ 



■r/ + (A^-i)2' 



(53) 



The density weighted two-vortex velocity field can then be 
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FIG. 3. Incompressible kinetic energy spectra for a single vortex 
(chain line), a vortex dipole (solid line), and a vortex pair (dashed 
line). The vortex dipole and pair are both shown for vortex separation 
d = 20^, and the wavenumber k^i = 2n/d is shown as a vertical 
dashed line. The cross-over scale Z:^ = 1 is given by the solid vertical 
line. 



written as 

/pi(r)p2(r) I I 

V = — [vi(r) + V2(r)] = Vpi(r)vi(r)+ VP2(r)v2(r) 



The final teim is 



Knir) 



Pi (r)p2(r) 



+^:i2(r). 



vi(r) 



(54) 



no 



P2(r) 



+V2(r) 



1 



no 



Pi(r) 



(55) 



which is only significant when considering the velocity field 
of one vortex in the close vicinity of the other vortex core. K12 
is therefore a negligible correction to the spectrum for widely 
separated vortices. This approximation may be trivially gen- 
eralized to arbitrary numbers of vortices provided the core s do 
not overlap appreciably. This approximation is central to our 
treatment, as it allows recovery of familiar point-like vortex 
physics in the infrared regime I. Strict validity is lim- 

ited to the regime where the intervortex spacing d is bounded 
below by ~ = In^, the scale at which the core structure 
becomes evident (see Fig.lTJ. 



d'r, 



We require the transform 



-ikr 



[ VP2(r)v2(r)]^. 



(56) 



where the superscript d denotes the case of a vortex dipole, 
and subscript j - x,y indicates the x and y components of the 
density-weighted velocity fields of each vortex. As above, the 
subscripts 1 and 2 denote vortices 1 and 2, and the superscript 
/ denotes that it is the incompressible, or divergence-free por- 
tion of the density-weighted field that is of interest here. To 
account for the opposite signs of circulation for the two vor- 
tices, without loss of generality we choose vortex 1 as posi- 
tively charged, located at ro - {d/2)x, so that Vi(r) = v(r-ro) 



where v(r) is the central vortex velocity field (|24] |. Vortex 2 
has velocity field V2(r) = -v(r 4- ro). We then have 



i- J ^2^e-'"c--°'[^v(r)j 



For a vortex dipole we can then wiite 

rf(k) = T/k)e"''' '"'= - :r;(k)e* ''° 

where f is the spectrum of a single vortex. Using 
the fact that 'Fyik) = -T^ik), we find for a vortex dipole 

E'lik) = i^^Fdkm - Mm- 



(57) 

(58) 
, and 

(59) 



The spectrum of a pair of vortices of the same circulation sep- 
arated by d is calculated similarly to be 



EP{k) = 2Q^Va(^^)(1 + Ukd)). 



(60) 

The spectra ( |59] l and (|60] | are shown in Figure |3] It is clear 
that for scales less than the vortex separation distance d there 
is interference in A:-space, leading to oscillations in the spec- 
trum. The difference between the dipole and pair is that the in- 
terference fringes are offset, and the infrared asymptotics are 
different, a feature we discuss further below. The spectrum of 
the vortex pair is clearly similar to that of the single vortex in 
the far infrared region, but the additional kinetic energy of the 
vortex pair state is observed throughout the spectrum. 



D. Kinetic energy spectrum of A'-vortex configurations 

Extending the above analysis, for a general system of 
singly quantized vortices with circulation signs Kp - +\ lo- 
cated at Tp the kinetic energy spectrum is given by 

N 

£f (^) = Q^VA(fc^) Yj '<pK,Mk\rp-r,\). (61) 

p=\,q=l 

We note the resemblance to point-vortex results which also 
have the Bessel function dependence ll29ll53ll . The function 
Ft^{k^) gives the incompressible limit for small k [Jo{kr) ap- 
proaches unity for small k and finite r], and introduces the 
physics of compressible superfluids for 1 < k^. 
We can write (1611 1 as 



where 



E^ik) = Nn^'F^{k^)GN{k) 



2 

GN{k) = 1 + '^P'^^MklTp - Tgl) 



(62) 



(63) 



p=l C]=p+l 



is a purely configurational function involving summation of 
M = A^(A^ - l)/2 distinct intervortex distances. This function 
has the limits 



lim G^ik) - 1 , 



lim 



(64) 



(65) 



p=l q=p+l 
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where the total dimensionless circulation is defined by 



IV. KOLMOGOROV SPECTRUM 



(66) 



for any contour C enclosing all positive and - N - N+ 
negative vortices. 

We then find that the onstrophy for the A^-vortex system is 



= Hm k^Efik) = A^A^Q = A^^,, (67) 



and consequently 



itf»i 



(68) 



This is one of our central results: the ultraviolet regime A:^ » 1 
has a universal asymptotic form that is independent of the vor- 
tex configuration, and that resembles the ultraviolet spectrum 
of classical 2D turbulence that is identified with a direct cas- 
cade of enstrophy. If we try to evaluate the classical definition 
of enstrophy, the result is singular, yet the onstrophy definition 
( |42] | gives a well-defined additive quantity that is singularity 
free and depends only on the total number of vortices in the 
system. 



E. Infrared behavior 



When r the far infrared limit 



gives 



E^{k)\ 



■3r2 



(69) 



This configuration-independent A: ' power law arises from the 
far-field velocity distribution of a collection of point vortices, 
which becomes equivalent to that of a single vortex of charge 
r at sufficiently large scales. 

When r s we use the small-argument expansion Jq{z) - 
1 — lA, and the asymptotic form (|49] l, to find the «: 1 
behavior determined by the configurational information con- 



tained in the intervortex distances |r,, 



(70) 



p=\ q=p+\ 



The simplest case involves a single vortex dipole and has only 
one length scale, namely the vortex separation, and the low-fc 
form E'jik) ^ Q^^d-k/2, as shown in Fig.|3] In general, when 
r = the infrared region of the spectrum is sensitive to the 
vortex configuration, but approaches a power-law for low-^ 
that has a configuration-independent exponent. The linear de- 
cay of kinetic energy as A: — > stems from the cancellation of 
the far-field velocity profiles for length scales greatly exceed- 
ing the largest intervortex separation in any neutral configura- 
tion of vortices. 



In the previous section we obtained an explicit expression 
for the incompressible kinetic energy spectrum that in- 
corporates the compressible nature of individual vortex cores 
through the function F^ik^) (derived via an ansatz for the vor- 
tex core profile), which captures the essential physics of the 
corresponding exact solution F{k^) defined in ( |36] |. For both 
functions, point-vortex physics is recovered at large length- 
scales (k^ «: 1). If the dynamical evolution is such that an 
inertial range associated with an inverse energy cascade devel- 
ops, we should expect a Kolmogorov power law Ei{k) oc k^^^^ 
over the inertial range. It is clear from the form of (l62i that 
this law can only depend on the spatial configuration of the 
vortices. We now seek to understand the simplest situations 
that may show evidence for the existence of such an inertial 
range. We consider forcing occurring via vortex and energy 
injection at a forcing scale kf ~ and describe vortex con- 
figurations that do and do not lead to a Kolmogorov law for 
k<^-^. 

We now assume an idealized case in which the infrared 
spectrum is continuous with the universal k^^ law of the ultra- 
violet spectrum at the scale k^ ^ 1. This constraint imposes a 
strong restriction on the infrared spectrum, completely deter- 
mining its form in the case that it satisfies a power law. In this 
respect the universal ultraviolet region has significant physical 
consequences. The power-law approximation to the universal 
ultraviolet region based on (l67T i has the form 



(71) 



The number of vortices determines the A^- vortex onstrophy 
( |67] l. from which the power-law approximation to the ultravi- 
olet energy spectrum is completely determined. This power 
law is a very good approximation, as will be seen by sam- 
pling different vortex configurations below. The infra-red or 
configurational regime is then given by the A:^ <K 1 regime of 
FAik^y. 



GN(k). 



(72) 



At this point we consider the consequences of assuming that a 
turbulent system will have a k^^^^ law in the configurational 
regime, and that this power law is continuous with (iTll at 
k^ - \. We suppose that E^^(k) oc k^^^^. Continuity at 

k - 1/^ then requires £^(1/^) = E^^(l/^), and gives the 
infrared spectrum 



A- 



/3 7,-5/3 



(73) 



Thus the constraint that the universal regime is continuous at 
the cross-over scale A: = completely constrains the form of 
the configurational spectrum. Physically this may correspond 
to an inertial range that extends upwards from the smallest 
scale of the infrared region given forcing at a wavenumber 
kf ~ Note that this expression ( |73] l has no reference to 
the signs of the vortex circulations, the degree of circulation 
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polarization, or vortex clustering. By assuming continuity at 
- I with a UV spectrum that has a universal A^- vortex form, 
we have implicitly assumed that all vortices are involved in 
determining the spectrum of the inertial range. 

We might expect that this will give a very good description 
for a completely polarized system exhibiting fully developed 
turbulence. When there is clustering in mixtures of different 
sign vortices, the spectrum may well still approach a Kol- 
mogorov law, but there is no reason to expect that it will cross 
over so smoothly. We return to this problem when we compare 
our analysis with numerical simulations of forced turbulence 
in Sec.lVCl 

It is useful at this point to give a simplified reiteration of 
Novikov's argument for the power law for the vortex distri- 
bution being -1/3 in Kolmogorov turbulence llssll . To obtain 
power-law behavior we must consider the spectrum for a vor- 
tex distribution involving many length scales. For simplic- 
ity we assume all vortices have the same sign of circulation, 
Kp = K. The configuration function ( |63l l has M - N(N -1)12 
terms in the summation, and can be written as 



GN{k) = 1 + "I X •'^^^'P^ 



(74) 



in terms of an average over distinct vortex separations s,. We 
introduce the intervortex distance distribution P{s) such that 
P{s)ds is the fraction of intervortex distances in the range 
[s, s + ds). In the continuum limit 



GN{k) oc P(s)Mks)ds. 



(75) 



We seek a distribution P{s) that will generate a Kolmogorov 
law from the A^- vortex spectrum (l62t for scales larger than the 
vortex core, k <s: ^ in the large-A^ regime. We then find 
from (|22]i that 



1 r°° 

- J P(s)Jo(ks)ds 



P{s)Jo(ks)ds. 



(76) 
(77) 



The scale invariance of turbulence naturally leads to the 
assumption that the intervortex separation distribution is a 
power-law P{s) ~ s^" over the scale range of interest. The 
requirement of a power law in the kinetic energy then gives 
the scaling relation 

1 1 r°° 
1 r°° 

k^-" Jo 



'Jo{T)dT. 



(78) 



The integral is convergent for -1/2 < a < I, allowing 

1 < j6 < 5/2. (79) 

In particular, the universal Kolmogorov law /3 - 5/3 occurs 
for 



P(s) 



-1/3 



(80) 



as obtained by Novikov 115311 for point-vortices. We will test 
this scaling argument for the exponent in numerical sampling 
of localized vortex configurations in the following sections. 
Testing if this vortex separation power-law holds in simula- 
tions and experiments may give a quantitative measure of fully 
developed 2D turbulence in a compressible superfluid, and a 
way to identify the inertial range as the scale range over which 
this power-law can be identified. 

In 2D classical turbulence a k^^ region of the kinetic energy 
spectrum is often associated with a direct enstrophy cascade. 
We note that this exponent jS = 3 is ruled out by (|79]l. Hence, 
within this continuum analysis the k^^ power-law spectrum 
cannot occur in the configurational region for 2D quantum 
turbulence, as long as the vortex distribution follows a sim- 
ple power law. This result suggests that if a direct enstro- 
phy cascade were to occur in the configurational region of 
the spectrum, a different type of vortex distribution would be 
necessary. This makes intuitive sense, since direct enstrophy 
cascades may be associated with the stretching of patches of 
vorticity in the 2D plane. Furthermore, direct enstrophy cas- 
cades and energy spectra proportional to k^^ have been noted 
in simulations of superfluid helium thin films ifsil Issll . Nev- 
ertheless, since we have also shown that the k^ > I range is 
determined entirely by the core structure, this gives a strong 
indication that a direct enstrophy cascade cannot occur in this 
region in compressible 2DQT. 

It is clear that configurations containing one or a few char- 
acteristic length scales, such as a vortex dipole or a vortex 
lattice, cannot lead to a power law spectrum for E^^f^ik). In the 
case of a vortex lattice the intervortex distance distribution has 
many discrete peaks ll36ll . The vortex dipole and a vortex pair 
each have a single length scale and this leads to characteristic 
interference fringes in the energy spectrum seen in Fig. |3] 



A. Sampling spatial vortex distributions 

We now test our analysis of the spectrum by numerically 
sampling several vortex distributions {rp,/<-p)^ j and evaluat- 
ing ( l62b . A straightforward test of the statistical argument 
for the Kolmogorov power law to occur involves sampling 
the power law dSOl l and evaluating (|74] i. Indeed, it is easily 
verified that this generates a k^^^^ spectrum in the configura- 
tional region. However, the connection of such a sampling to 
particular spatial vortex distributions is not clear, and in fact 
the mapping is not unique. To make this connection concrete 
we require a way of sampling finite, localized, spatial vortex 
distributions for which the vortex separations are power-law 
distributed. 

In an ideal, infinitely extended vortex configuration exhibit- 
ing the power law (ISOl l. the system is translationally invariant 
and the coordinate origin can be placed at any particular vor- 
tex, yielding the same power law for the radial distribution of 
vortices from the origin. In a finite system the scale invariance 
can only persist up to scales of order the largest vortex sepa- 
ration. Furthermore, the vortices must be separated by a min- 
imum distance to satisfy the assumptions used in deriving the 
energy spectrum from the point-vortex model. In practice it 
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is necessary to use a self-consistent sampling scheme in order 
to generate the correct power-law distributions for localized 
finite configurations. 

Our sampling scheme for a configuration of vortices is: 

1 . Sample the radial distance Vp of each vortex from the 
coordinate origin according to a power law probability 
distribution oc r^f . The exponent a is distinct from a 
due to the finite system size and localization of the dis- 
tribution. 

2. Assign each vortex a randomly chosen, uniformly dis- 
tributed angle 9p € [0, 2n). The cartesian coordinates 



for vortex p are then 



{xp,yp) = rp(cos0,„sin0;,). 



(81) 



3. For a given A^, compute the kinetic energy spectrum, 
averaging over n 5 = 100 samples ofvortex position data 
found via the foregoing routine. Iterate until the spectral 
power law of interest is found. 

We will sample a power-law distribution of vortex distances 
from the origin. In practice there is a lower (r^m ~ ^) and 
upper (rniin ~ R = system size) cutoff for power-law scaling. 
We thus wish to sample the distribution 



Pair) 



1 



that is normalized on the interval r^in ^ 
r values from this distribution usin: 
X e [0, 1 ) via the transformation ||5( 



(82) 

< r < Tmax- We sample 
uniform random variates 



rmax^ ' min 



n'(l-^)l' 



(83) 



Note that when a < 1, as is always the case in this work, 
some kind of ultraviolet cutoff is required for the distribution 
to be normalizable. Here we have made a choice that gives 
a power-law distribution over a well-defined scale range (See 
e.g. Ref. lf56ll for other common choices). 



B. Classical velocity distribution of a large cluster 

A fully polarized configuration of vortices with a given ra- 
dial power law forms a quantum analogue of the coherent 
vortices of forced 2D turbulence in classical fluids. As an ar- 
rangement of many vortices, the velocity distribution must ap- 
proach a classical limit, according to Bohr's correspondence 
principle, in much the same way that a rotating Abrikosov lat- 
tice generates a velocity field that approaches that of a rotating 
rigid body jsoll . In what follows we find the classical veloc- 
ity field, and identify the physical significance of the radial 
exponent a. 

To determine the velocity field we first compute the fraction 
of vortices enclosed by a circular contour around the origin. 



with radius r, for the distribution (|82| |. Taking r^in - ^, rmax = 
R, and considering scales ^ <s; r <sc /? we have 



where we used the fact that - s~"ds - (1 - 

a)-\R^-" - ^'-*) ^ R^-"/(l - a) normalizes the distribution 
up to the largest scale R. Considering the average azimuthal 
velocity component v^{r), the circulation is 



■ dl - —n - 

m 



'r)2nr 



(85) 



where « - f,-N is the number of vortices enclosed by the con- 
tour of radius r. Using ( [84l ), we obtain the velocity profile 



(0- 



cN 



{Ri^y-%ri^r 



(86) 



This inertial cluster has a power-law velocity profile deter- 
mined by the specific radial exponent a. In contrast, the veloc- 
ity profile of an Abrikosov vortex lattice rotating at frequency 
0) approaches that of a rigid body M^(r) - cor, and a single 
quantum vortex has profile M^(r) - H/mr. The inertial cluster 
velocity profile (|86] | is compared with sampled distributions 
(see below) in Figure |4] 



C. Scale expansion of a large cluster 

We now illustrate the role of the radial power law expo- 
nent a and the classical velocity distribution by sampling a 
large cluster We vary our choice of scale range for the vor- 
tices (rtnax in Eq. [82] |. and investigate how the range of k^^^^ 
changes. 

A characteristic wave number measuring cluster size for a 
sample involving iic vortices in a given cluster is given by 



Where r 



k = Inl? (87) 
In the figure we plot the scale k to 



give an indication of the range of k^^^^ scaling. We also plot 
^max corresponding to the largest vortex separation scale in the 
system 



^max = 27r/max|r„ - 



(88) 



For k < ^max the velocity field approaches that of a single 
vortex. 

In Figure |4] vortex distributions are sampled for = 100 
vortices of the same sign using the sampling scheme dSTT l. 
(l83T l, for r,nin = f and different values of rmax- Individual sam- 
ples are shown to indicate the spread of vortices, and the ve- 
locity profiles and kinetic energy spectra are computed by av- 
eraging over n, = 100 samples. The mean azimuthal velocity 
compares well with Equation ( [86] l, showing that the specific 
radial exponent a in (|82] i also determines the power law of 
the azimuthal velocity field, as seen in ( |86] |. For scales larger 
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FIG. 4. Locations of N = 100 vortices for particular samples (left), with velocity distributions [middle, with the analytical result ( I86t (dashed 
line)], and kinetic energy spectra (right, blue curves), averaged over n, = 100 samples distributed according to ( I81t and ( I82t . with = <f, 
and (a) r„y^^ = 60n^, a = 0.4, (b) r,nax = 400;r<f, a = 0.65, and (c) r^ax = 4000;r<f, a = 0.83. The increasing scale range of power-law behavior 
gives an inertial range in the corresponding kinetic energy spectrum of ~ 1, 1.5, and 2 decades respectively. 



than Tniax, ^vCO returns to the r"' scaling for a charge vortex 
[clearly seen in Fig.|4](fl)]. = 100 vortices distributed up to 
'"max - 6Qn^ gives an inertial range in the kinetic energy spec- 
trum of approximately one decade, for a - 0.4 (as the number 
of vortices in a given scale range increases, a ^ 1 /3). In- 
creasing the upper scale cutoff to r,nax - 4007r^, 40007r^ gives 
a = 0.65,0.83, with 1.5 and 2 decades of inertial range re- 
spectively. Thus expanding the scale range of power-law be- 
havior expands the inertial range, but requires the azimuthal 
velocity profile to steepen. In the UV region of the spectrum 
the law always holds [Eq. (fTTTil, while the inertial range 
holds for ^ < A: < ^ and single vortex behavior is apparent 

for k < ^max- 



D. Clustering in a neutral distribution 

We now consider the role of increasing clustering in ex- 
panding the inertial range of the kinetic energy spectrum. 



First, in «j = 100 samples, we distribute = = 100 
vortices over a uniform periodic domain; one such sample 
is shown in the inset of Figure |5] (a). The corresponding ki- 
netic energy spectrum shows the correct UV-region spectrum 
given by Eq. (fTTT i. and also approaches the Ei(k) cc k form for 
k < kmax, given by Eq. (iTOl i. For fc^ax ^ k < the spec- 
trum is less steep than k^^^^ and the system lacks an inertial 
range. In Figure |5](b) the vortices are sampled as 40 clusters 
of «(, = 5 vortices of the same sign according to (ISTt and (l83T l, 
with a - 0.8, and r,nin = ^, rmax - lOn^. The 10 +ve and 10 
-ve cluster centers are uniformly distributed as in Fig. 0(a), 
as seen in the sample (inset). This distribution yields a k^^^^ 
power-law kinetic energy spectrum over ~ 1 decade of wave 
numbers. By further expanding the scale of clustering while 
reducing the number of clusters to preserve and A^_, we 
find the inertial range can be extended. In Figure |5](c) sam- 
ples consists of lie - 20 vortices in each cluster, distributed 
between r^i„ - and r^ax - IOOtt^, with a = 0.75, and 
giving ~ 2 decades of inertial range. We note that, compared 




FIG. 5. Kinetic energy spectra for A'+ = 100, and N- = 100 vortices, 
averaged over n., = 100 samples (blue curves). The vortices are dis- 
tributed (a) uniformly over the (600^)- square domain (inset), and 
clusters are sampled according to ( 181) and ( 182) , with (b) = 5 vor- 
tices in each cluster, r^in = r,i,ax = Wn^, and (c) = 20 vortices 
in each cluster and = 100;r^, with each cluster center uniformly 
distributed over the periodic domain as in (a). The power-law sam- 
pling in figures (b) and (c) requires a radial exponent of o' = 0.8, 
and 0.75, and gives a Kolmogorov k^^^^ region in the corresponding 
kinetic energy spectrum of ~ 1 and 2 decades respectively. 



with Fig.m k, shown in Fig.|5](b) and (c) does not correspond 
so well with the lower bound on the inertial range, presum- 
ably because of the significant space between clusters in the 
neutral system. 



V. FORCED 2D QUANTUM TURBULENCE 
A. Scenario of forced turbulence 
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FIG. 6. Illustration of an inertial range (non-shaded region) for 
the incompressible portion of kinetic energy in forced compressible 
2DQT. The E(k) oc k^^ region arises from the structure of the vortex 
core and thus is not a signature of vortex configurations and vortex 
turbulence. This ultraviolet region can thus not support energy cas- 
cades, nor does this region correspond to enstrophy cascades. Net 
energy injected at kf ~ in the form of vortices can only move to- 
wards the infrared. The Kolmogorov law E(k) oc k^^^^ occurs in the 
inertial range of fully developed turbulence. The far-infrared region 
is given by E(k) oc A; for a system with no net vorticity, and is evident 
for k <K L"' where L is the largest intervortex distance. For forcing at 
smaller wavenumbers, the spectrum may be more complex, possibly 
involving other forms of energy and enstrophy flux. 



in a compressible 2D superfluid 

The canonical model of 2D classical turbulence consists of 
a velocity field described by the 2D Navier-Stokes equation 

^ -H(v V)v = --Vp-H vV^v-^v + fv. (89) 
of p 

The density p of the incompressible fluid is held constant 
by the pressure field p, v is the kinematic viscosity, fy is a 
forcing term and A represents linear frictional damping aris- 
ing from irreducible 3D aspects of the system in which the 
2D flow resides. If the fluid is subjected to suitable forcing 
it will develop an inverse energy cascade and a direct enstro- 
phy cascade, with associated k^^^^ and k^^ power laws respec- 
tively lfl3ll . An inverse energy cascade induced by small-scale 
forcing can be steady because the -A\ term damps energy at 
large length scales Jgl]. For a homogeneous compressible su- 
perfluid subject to forcing from an external potential, (ITZt can 
be written as 

^H-(vV)v = --Vp-i-v„V-v-Hfv, (90) 
ot m 

where the forcing fy = -Vy(r, ?)/'« is assumed to be spa- 
tially localized. The lack of a -A\ frictional damping term 
means, in the classical case, that if an inverse energy cascade 
develops as a result of steady forcing, it is not expected to be 
stationary. In the superfluid case, the compressibility of the 
fluid allows vortex-antivortex annihilation, which couples en- 
ergy into the sound field. This interaction between the sound 
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and vorticity fields renders the calculation of energy fluxes 
in compressible superfluid systems particularly difficult and 
somewhat ambiguous llsoll . 

As we have shown above, in contrast with the classical 
Kraichnan scenario utilizing a 2D Navier-Stokes analysis, a 
k^^ spectrum for A: > ^ ' for a 2D quantum fluid is not caused 
by a direct enstrophy cascade but is rather a consequence of 
the vortex core structure, and thus should not be interpreted 
in terms of vortex configuration dynamics (note that vortex 
core shape excitations can be neglected, since they constitute 
a component of the sound field). Given forcing at a wavenum- 
ber kf ~ and minimal vortex-antivortex annihilation, the 
incompressible kinetic energy can only move toward the in- 
frared. This scenario is shown schematically in Figure |6] 

It has been shown that dipole recombination provides a 
route for a direct energy cascade to develop in 2D GPE dy- 
namics isoll . This mechanism provides a means for opposite- 
sign vortices to approach zero distance, coupling vortex en- 
ergy to the sound field during vortex annihilation. However, if 
the forcing leads to significant clustering of like-sign vortices 
faster than recombination occurs, or prior to recombination 
occurring, dipole-decay will be strongly inhibited. This sug- 
gests that under the right conditions of forcing an inverse en- 
ergy cascade can become the dominant mechanism of energy 
transport between distinct length scales. 



B. Kolmogorov constant and clustered fraction 

By making use of the universal onstrophy and the condition 
of continuity atk^ ^ 1 we have found that the k^^^^ power-law 
given by (l73T l describes the spectrum of numerically sampled 
vortex configurations that exhibit a s^^^^ power law for the 
vortex separation data. While individual spectra and configu- 
rations do not give information about dynamics, in particular, 
the direction in A:-space of any energy cascades, the power 
law suggests the existence of an inertial range comprised of 
vortices. In a cascade, such a configuration will transfer in- 
compressible energy between scales while conserving energy. 
Assuming the infrared portion of our double-power-law anal- 
ysis, namely ( |73] l, will also describe such a cascade, we can 
cast it as a statement about the Kolmogorov constant in terms 
of the one-vortex onstrophy and the slope of the radial wave- 
function at the vortex core. 

To write ( |73] l in standard form, we introduce the unique A^- 
vortex quantity with dimensions energy/mass/time that can be 
constructed from Q^^, m, and ti/Q.^^: 



We then find 



mh 



(91) 



(92) 



where the remaining quantities have been absorbed into the 
dimensionless Kolmogorov constant: 



1/3 



ju = tP'/m^^ in the homogeneous system, and the bar notation 
distinguishes the quantum system. In the dilute Bose gas the 
2D interaction parameter is fi/ no - g2 = Antra /ml, where 
/ is the characteristic thickness of the three-dimensional sys- 
tem I30L|57|1 . In terms of this length we find 
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(94) 



We emphasize that the physical input needed to arrive at 
this form of the Kolmogorov constant is (i) accounting for 
the structure of a compressible quantum vortex in determin- 
ing the ultraviolet spectrum, and (ii) imposing continuity of 
the ultraviolet spectrum at A:^ = 1 to a Kolmogorov power- 
law in the infrared. In classical turbulence, €20 - 1 IdSII . To 
give an example of how €20 may be evaluated for a compress- 
ible superfluid exhibiting 2DQT, we consider a '^^Rb BEC that 
is homogeneous in the x - y plane, and that is harmonically 
trapped in the z-dimension with trap frequency u- - 27rx5000 
Hz. We use m - 1.44 x lO^^^kg, a = 5.8nm, for which 
I - yjlntilmuz - 0.38 //m, and g2 - 0.197tf/m. For these 
values, €20 - 0.212. Note also that by defining the configura- 
tional rate constant ( |9TT i. we have confined this discussion to a 
scale invariant distribution involving vortices. This expres- 
sion suggests that such a configuration can support an inverse 
energy cascade at the rate e^v. 

The foregoing discussion involves an ideal distribution of 
vortices configured with the ff = 1 /3 power law. It is clear 
from Fig.|5]that the vortices do not have to all have the same 
sign of circulation, but they must be configured into clusters 
of vortices with the same sign. The universality of €20 in 
classical turbulence in incompressible fluids leads us to postu- 
late that the condition that all vortices are power-law clus- 
tered can be further relaxed. For fully developed quantum 
turbulence involving A^ vortices, we interpret A^ as the par- 
ticipation number, representing the number of vortices in a 
scale-free turbulent configuration, which in the case of a fully 
polarized cluster is maximal. Imperfect clustering involves 
fewer vortices in power-law cluster configurations, and an ef- 
fective participation number that is the number of clustered 
vortices A^^, < A^, namely, the number with nearest-neighbors 
of the same sign. Making the replacement A^ — > A^^. in (73[ , 
we propose the ansatz spectrum 



(95) 



(93) 



as a more general definition for systems that have incomplete 
clustering in the inertial scale range. We test this hypothesis in 
the next section in dynamical simulations of the forced dOPE. 

It is important to note that the condition of continuity at 
k^ - I is no longer exactly met, since only the value A^^. 5 A^ 
will produce an infrared spectrum with k'^^^ that is continuous 
with the ultraviolet power law approximation at k = 1 /^. We 
also note that a more general measure of clustering, namely 
the polarization index (P) was introduced in Ref. 15^ mea- 
suring the degree and type of spatial clustering of like-sign 
vortices in 3DCT. The Kolmogorov k^^^^ spectrum was found 
to correspond to the partially polarized value P = 1/3, while 
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FIG. 7. Time evolution of grid turbulence in damped GPE. Left: particle density (rescaled to the peak density). Center: vortices colored by 
charge, with total (A') and clustered (A'^.) numbers of vortices. The tield of view is (1024^)^. Right: Incompressible energy spectra (circles), 
with the Kolmogorov ansatz [red line, Eq. l|95ll1, the ansatz for a polarized cluster of A' vortices [dashed line, Eq. ( I73H , and the universal k^^ 
region [blue line, Eq. M\\\. 



Other scaling laws yield differing polarizations. The clustered 
fraction used in the present work is a simpler (global) mea- 
sure of polarization as it does not contain information about 
the spatial distribution of vortices, but it is only relevant for 
2DQT. 

We can also write down an expression for €20 for a gen- 
eral energy spectrum that may be computed numerically from 
simulation data, by making use of the ansatz ( |95] l. This is 
equivalent to using (l92T i with — > A^^., from which we can 
define the function 



1/3 



(96) 



In a region where the spectrum is approximately k ^1^, ClAk) 



will be approximately constant, and it may be compared with 
the prediction (l93T l. We test this numerically in the next sec- 
tion. 

We note that in 3DQT an energy bottleneck has been pre- 
dicted for the direct energy cascade llsoi I60I1 . and also ob- 
served in GPE simulations llilll . It occurs due a mismatch 
between the rates of energy transport at large length scales 
(hydrodynamic regime) and small length scales (Kelvin wave 
cascade). The mismatch causes energy to pile up at the length 
scale where the two cascades meet. This raises the possibility 
of a bottleneck in 2D, although Kelvin waves are disallowed 
in 2DQT, so this particular mechanism would not be relevant. 
However, for a given forcing mechanism, it is possible that the 
rate of transporting energy to large length scales may not be 
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FIG. 8. Plot of C2p(k) [see Eq. | |96H computed from the grid turbu- 
lence kinetic energy spectra of Figure |7] The horizontal dashed line 
gives CiD from Eq. 04\ . 



high enough to remove all of the vortex energy introduced at 
the forcing scale. Thus a bottleneck could still occur, and our 
assumption of continuity of the spectrum at ^ 1 may not 
hold in general. We return to this question in the next subsec- 
tion, where we find some indication of an energy bottleneck 
at the forcing scale in numerical simulations. 



stirring beams. A small amount of initial noise is added to the 
wave function to break the reflection symmetries of the sys- 
tem. We thus evolve the system according to ([TJ for the same 
potential, but with the Galilean transformed nonlinear oper- 
ator X. ^ -C + ift^odx- During evolution the dimensionless 
damping rate is set to y = 0.003. 

The time evolution of the system is shown in Figure Q at 
three times, at approximately (1 /3, 2/3, l)L/vo, so as to avoid 
any periodic flow effects in the jc-direction. The four obstacle 
beams generate many vortices (up to N ~ 10^), and significant 
clustering (Nc/N > 0.6). The incompressible kinetic energy 
spectrum shows a wide region that is well described by the 
k^^^^ form of Equation ( |95] l. For later times [Figure |7](b), (c)] 
the spectrum shows a significant pile up around the forcing 
scale kf ~ suggesting a mismatch between the rates of 
injection and transport of incompressible kinetic energy. 

In Figure |8]we compare the function C^jjik) [Eq. (|96]l1, as 
numerically computed from our simulation data, with the ana- 
lytical prediction of the Kolmogorov constant €20 [Eq- (I94li1. 
The region of k^^^^ appears as a broad flat region that is in 
close agreement with €20 - 0.212 pertaining to our simula- 
tion parameters. 



VI. CONCLUSIONS 



C. Damped Gross-Pitaevskii Dynamics 

We now consider a simulation of the forced dGPE that gen- 
erates significant clustering of vortices of the same sign. The 
system consists of a homogeneous superfluid with periodic 
boundary conditions, stirred by dragging four Gaussian ob- 
stacle beams through it at a constant speed |32|. thus model- 
ing grid turbulence in a BEC. When an obstacle is dragged 
through a superfluid sufficiently rapidly, superfluidity cannot 
be maintained. For slowly moving obstacles, the superfluid 
will adapt to the forcing, and vortices are not formed. Above a 
critical velocity v^ ll62ll vortex dipoles are periodically formed 
in the wake of the obstacle, injecting linear momentum into 
the superfluid. Sufficiently rapid motion (v » v^) causes 
many vortices to be nucleated behind the obstacles in a chaotic 
fashion ll63ll . involving clustering of like sign vortices. Our 
choice of obstacle speed puts the system dynamics in the lat- 
ter category. 

We work in units of jj, ^, and ^/c for energy, length, and 
time, respectively. In these units the specific parameters we 
choose (see the previous subsection) are g2 = 0.197//^-, cor- 
responding to homogeneous density no = f^/gi - 5.26^"^, 
and Niot = 5.5 x 10*" particles in a homogeneous 2D system 
of side length L = 1024^. The Gaussian potentials each have 
fixed 1/e^ width of wo - VSf, and height Vo - 100//, and are 
initiafly located at x = -L/2 + 8^,y - ±L/8, +3L/8. Numer- 
ically, we proceed by first finding a ground state of dl}, for a 
homogeneous system with periodic boundary conditions, sub- 
ject to the localized obstacle beams. We then transform into 
a frame translating at vq = 0.8c, and maintain the obstacle 
locations relative to this frame, creating a dragging grid of 



To summarize, we have investigated relationships between 
the concepts of 2D turbulence in classical fluids and the 
emerging topic of 2D quantum turbulence of vortices, specif- 
ically as it relates to Bose-Einstein condensates. We estab- 
lished a link between the hydrodynamic limit of the damped 
GPE and the Navier-Stokes equations, providing an estimate 
of a quantum Reynolds number for superfluid flows in BECs. 
We have given a theoretical treatment of the incompressible 
kinetic energy spectrum that explicitly incorporates the vor- 
tex core structure in a compressible superfluid. The incom- 
pressible kinetic energy spectrum for a compressible super- 
fluid is deconstructed in terms of single-vortex contributions 
determining a unique ultraviolet power-law where the energy 
spectrum scales as k^^, and a contribution that depends on 
the configuration of vortices within the fluid that determines 
the infrared region of the spectrum. For the configurational 
regime we find: 

1 . The spectrum only depends on the distribution of vortex 
separations and the sign of the circulation of each quan- 
tum vortex. If the distribution of vortex separation s for 
a system of vortices of the same sign is a power-law 
oc 5 " with exponent a = 1/3, the kinetic energy spec- 
trum will take the universal Kolmogorov form oc k^^^^, 
as shown for point vortices ll53ll . Localized clusters of 

vortices of the same circulation with this power law 
distribution can be constructed by sampling using a spe- 
cific radial exponent a that depends on the number of 
vortices and the scale range over which they are dis- 
tributed. 

2. The azimuthal velocity field of a large cluster is deter- 
mined by a. By inflating the scale range of a cluster we 
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find that a increases, the velocity field is steepened and 
the inertial range expands to larger scales. In a neutral 
system the inertial range can be extended by increasing 
the size of clusters while decreasing their number 

3. The universal form of the UV region of the kinetic en- 
ergy spectrum imposes a strong constraint. If the Kol- 
mogorov power law occurs in the infrared region, then 
the postulate of continuity between the infrared and ul- 
traviolet regions completely determines the spectrum 
when the ultraviolet and infrared regions are approxi- 
mated as power laws. Physically, this corresponds to 
the inertial range extending down to the smallest con- 
figurational scale of the system ~ ^. We note that the 
postulate of continuity may not be relevant for all sys- 
tems or forcing mechanisms. 

4. We infer an analytical value for the Kolmogorov con- 
stant [Eq. ( l93T l1 under the conditions of spectral conti- 
nuity at the cross-over scale for a system of vortices of 
the same sign. To assess the validity of this inference 
for dynamical situations we compare our analytical re- 
sults with spectra from a numerical simulation of the 
forced dGPE for the specific case of a dragging a grid of 
obstacles through an otherwise homogeneous BEC. We 
find reasonable agreement provided we introduce the 
concept of a clustered fraction Nc/N < 1, which is the 
fraction of vortices that have same-sign nearest neigh- 
bors. This measure discounts all vortex dipoles from the 
configurational analysis. We then observe good agree- 
ment between our Kolmogorov ansatz [Eq. (I95ll1. and 
the spectrum calculated from the dGPE data. We also 
find that the predicted value of the Kolmogorov con- 
stant is in close agreement with the numerical simula- 
tions [Figure [8l. 

We note that while our analysis indicates that vortex positions 
and circulations are enough to determine an approximate in- 
compressible kinetic energy spectrum, the reverse is not nec- 
essarily true: a Kolmogorov spectrum does not carry infor- 
mation about specific vortex distributions. Nevertheless, our 
analysis does indicate that the number of vortices in a quantum 
fluid can in principle be directly determined from the ultravi- 
olet energy spectrum. Moreover, the concept of a cascade in 
turbulence implies system dynamics and energy transport, yet 
aside from our numerical simulation example, our analytical 
approach is an instantaneous measure. Importantly, one must 
determine means of characterizing vortex motion and relate 



such dynamics to the cascade concept. 

The field of 2D quantum vortex turbulence is relatively new, 
compared with the much longer histories of 3D superfluid tur- 
bulence, 2D classical turbulence, and even dilute-gas Bose- 
Einstein condensation. Point-vortex models have been ex- 
tensively used in descriptions of superfluid dynamics as well 
as in 2D classical turbulence, although point-vortex distribu- 
tions can only serve as approximate models of real 2D classi- 
cal flows. Our approach merges concepts from each of the 
above subjects in order to develop a new understanding of 
2D quantum turbulence. By considering the compressibility 
of a dilute-gas BEC, we find an analytical expression for the 
ultraviolet incompressible kinetic energy spectrum and an A^- 
vortex equivalent of enstrophy in a quantum fluid, which we 
term the onstrophy. For the infrared region, point vortex mod- 
els are sufficient, and vortex configurations serve to identify 
spectra as summarized above. Taken together, the primary 
new outcome of our analysis is a link between vortex distri- 
butions, vortex core structure, and power-law spectra for 2D 
compressible quantum fluids. 

Future work on 2D quantum vortex turbulence will involve 
numerical simulations and comparisons with our analytical re- 
sults, extension of this analysis to confined and inhomoge- 
neous density distributions, characterization of vortex dynam- 
ics and the time-dependence of energy spectra particularly in 
relation to clustering f64|[, inverse-energy cascades, and the 
nonthermal fixed point [|3lll65ll . and investigation of connec- 
tions with weak- wave turbulence in BEC ifia fisj [l9i 166^6911 . 
We also believe that observing vortex distributions such as the 
ff = 1 /3 power-law for localized clusters may provide a new 
means of quantitatively characterizing 2D quantum vortex tur- 
bulence through direct experimental observations of vortex lo- 
cations in a forced 2D superfluid, and we are working towards 
realizing such experimental observations. 
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